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Abstract 

A  simplified  model,  to  that  produced  previously  by  the  authors,  for  the  galvanostatic  discharge  of  primary  alkaline  battery  cathodes  is 
presented.  Laplace  transform  and  perturbation  methods  are  employed  to  obtain  the  leading  order  spatial  and  temporal  behaviour  of  the  porous 
cathode  over  two  distinct  size  scales.  It  is  shown  that  for  a  wide  range  of  industrially  relevant  discharge  conditions  the  time  taken  for  KOH 
electrolyte  to  diffuse  into  a  porous  electrolytic  manganese  dioxide  particle  is  fast  compared  with  the  cathodic  discharge  time  and  that  ohmic 
losses  within  the  graphite  phase  of  the  cathode  can  be  considered  to  be  negligible.  Numerical  solution  of  the  simplified  model  equations  is 
discussed  and  the  results  are  validated  against  relevant  experimental  data. 

©  2005  Elsevier  B.V.  All  rights  reserved. 
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1.  Introduction 

Primary  alkaline  battery  cathodes  are  compacted  mix¬ 
tures  of  fine  particles  of  electrolytic  manganese  dioxide 
(EMD)  and  graphite  that  are  flooded  with  concentrated 
(7-10  M)  aqueous  potassium  hydroxide  (KOH).  Cathodes 
are  usually  constructed  with  relative  weight  percentages  of 
approximately  80-87%  EMD,  5-12%  graphite  and  6-12% 
KOH  [1].  These  cathodes  are  complex  porous  structures, 
and  their  discharge  behaviour  is  the  result  of  interconnected 
and  non-linear  physical,  chemical  and  electrochemical 
phenomena  that  occur  over  several  size  scales.  Development 
of  improved  battery  materials  and  systems  requires  a  sound 
understanding  of  these  processes. 

Mathematical  modelling  provides  a  unique  opportunity  to 
understand  the  processes  and  their  interactions  that  determine 
the  discharge  behaviour  of  primary  alkaline  battery  cathodes. 
Such  modelling  can  also  provide  direction  for  further  research 
and  explanations  for  experimental  observations. 
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In  a  series  of  papers,  Newman  and  co-workers  [2-8] 
developed  an  approach  to  modelling  porous  electrodes, 
referred  to  as  Macrohomogeneous  Porous  Electrode  Theory. 
This  approach  has  been  used  as  a  basis  for  a  number  of 
alkaline  battery  models,  including  those  of  Wruck  [9],  Chen 
and  Cheh  [10,11],  Podlaha  and  Cheh  [12,13],  Zhang  and 
Cheh  [14],  Farrell  et  al.  [15],  and  Kriegsmann  and  Cheh 
[16]. 

To  date,  the  systematic  analysis  of  previous  models  has 
been  largely  restricted  to  sensitivity  arguments  that  involve  a 
small  sub-set  of  parameters  within  the  model  system  in  ques¬ 
tion.  A  more  general  sensitivity  analysis,  however,  has  been 
undertaken  by  Zhang  and  Cheh  [17]  for  the  alkaline  battery 
model  originally  posed  by  Chen  and  Cheh  [  10] .  Furthermore, 
Farrell  and  Please  [  1  ]  applied  perturbation  methods  in  order 
to  analyse  a  set  of  equations,  previously  developed  by  Farrell 
et  al.  [15],  for  the  galvanostatic  discharge  of  a  single  porous 
particle  of  EMD.  Their  work  led  to  significant  simplifications 
and  a  greater  knowledge  of  the  behaviour  of  particle  dis¬ 
charge.  Nevertheless,  we  note  that  no  experimental  validation 
of  their  results  was  presented  because  individual  particle  dis¬ 
charge  at  the  size  scales  found  within  commercial  cathodes 
is  difficult. 
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Nomenclature 

a 

activity  of  a  dissolved  species  (mol  cm-3) 

A,B 

arbitrary  constants 

C 

concentration  (mol  cm-3) 

Qvin4+ 

Laplace  transform  of  CMn 4+ 

Ce.0,1... 

terms  of  the  asymptotic  expansion  of  Ce  in  (*5 

^eoo 

diffusion  coefficient  of  KOH  (cm2  s-  1 ) 

/ 

function  used  to  solve  auxiliary  Eq.  (40). 

F 

Faraday’s  constant  (C  mol-1) 

H 

height  of  cathode  (cm) 

i 

current  density  (A  cm-2) 

in 

transfer  current  density  from  oxide  to  solution 
phase  (A  cm-2) 

*0 

initial  equilibrium  exchange  current  density 
(A  cm-2) 

/ 

total  discharge  current  (A) 

k-2,  £3 

constants  associated  with  calculating  the  con¬ 
ductivity  of  EMD 

L 

the  Laplace  transform  function 

P 

Laplace  transform  variable 

r 

radial  coordinate  on  the  particle  scale  (cm) 

r0 

outer  radius  of  particle  (cm) 

R 

radial  coordinate  on  the  cathode  scale  (cm) 

^gas 

universal  gas  constant  (J  K- 1  mol- 1 ) 

Ri 

inner  radius  of  cathode  (cm) 

Ro 

outer  radius  of  cathode  (cm) 

t 

time  (s) 

t* 

dummy  variable  used  in  integrals 

rK+ 

transport  (transference)  number  for  K+ 

T 

temperature  (K) 

vD 

volume-average  velocity  of  electrolyte 
(cm  s-1) 

% 

partial  molar  volume  of  electrolyte 
(cm3  mol-1) 

Fii2o 

partial  molar  volume  of  H20  (cm3  mol-1) 

^Mn3+ 

partial  molar  volume  of  Mn3+  (cm3  mol-1) 

^Mn4+ 

partial  molar  volume  of  Mn4+  (cm3  mol-1) 

y 

radial  coordinate  on  the  crystal  scale  (cm) 

y0 

outer  radius  of  EMD  crystal  (cm) 

Greek 

a  1-17 

dimensionless  constants  associated  with  the 
dimensionless  equation  set 

£emd 

fraction  of  total  cathode  volume  that  is  EMD 

£s 

fraction  of  the  cathode  volume  that  is  void 

1 1 

local  overpotential  (V) 

Keoo 

conductivity  of  bulk  electrolyte  solution 
(S  cm-1) 

^•n 

the  positive  roots  of  Eq.  (28) 

conductivity  of  bulk  graphite  (S  cm-1) 

Superscripts 

0 

value  at  t  =  0  s 

Subscripts 

e  electrolyte  phase 

Mn4+  manganese  4+  ion 

p  manganese  oxide  particle  scale 


In  this  work,  we  begin  by  considering  the  full  cathodic 
discharge  model  proposed  by  Farrell  et  al.  [15],  which  was 
shown  to  compare  favourably  with  available  galvanostatic 
discharge  data.  Our  aim  is  to  apply  perturbation  methods  in 
order  to  analyse  the  key  physical,  chemical  and  electrochem¬ 
ical  processes  that  govern  the  behaviour  of  the  model  and, 
in  so  doing,  obtain  a  simplified  model  system  that  accounts 
for  these  key  processes.  This  work  extends  the  work  of 
Farrell  and  Please  [1]  into  the  cathodic  domain  and  yields  a 
simplified  model  with  solutions  that  can  be  directly  validated 
against  experimental  data.  Such  validation  is  presented 
here. 


2.  Model  development 

2.1.  Full  cathode  model 


To  model  effectively  the  processes  that  occur  in  the  multi¬ 
scale  environment  of  an  alkaline  battery  cathode,  Farrell  et  al. 
[15]  employ  a  simplified  description  of  the  cathode.  In  their 
model,  the  EMD  and  graphite  in  the  cathode  are  assumed  to  be 
well  mixed,  that  is,  each  porous  EMD  particle  is  surrounded 
by  a  continuous  graphite  phase  that  connects  the  porous 
EMD  particles  to  the  current-collector.  These  porous  EMD 
particles  consist  of  aggregates  of  non-porous  EMD  crystals. 
Both  crystals  and  particles  are  assumed  to  be  spherical.  In 
addition,  the  network  of  pores  within  the  EMD  particles  and 
between  the  EMD  and  graphite  particles  is  assumed  to  be 
flooded  with  concentrated  aqueous  KOH  electrolyte  solution. 
A  schematic  diagram  of  these  size  scales  is  given  in  Fig.  1. 
Furthermore,  the  Mn4+  ions  within  the  crystals  are  assumed  to 
only  undergo  the  first  electron  reduction  step  as  given  by  the 
reaction. 


Mn02  +  H20  +  e“  MnOOH  +  OH-  (1) 


The  governing  equations  and  associated  boundary  and  ini¬ 
tial  conditions  previously  introduced  by  Farrell  et  al.  [15]  to 
describe  the  galvanostatic  discharge  of  the  three-scale  struc¬ 
ture  depicted  in  Fig.  1  are  summarised  below.  The  equations 
are  given  in  dimensionless  form  below. 

The  crystal  scale: 


9^Mn4+  _  ( ~'2  ^Mn4+  ^ 

1  dt  ~  y2  dy  y  dy  )  ’ 


(2) 
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Scale  1 :  Cathode 
Typical  thickness:  1-5  mm 


Scale  2:  Porous 
EMD  particle 


C'e(p)  |  i  —  C, 
^(P)|r=l  = 

9Ce(p) 


dr 


®(P) 


r=0 


dr)(p) 


dr 


=  0, 

f= o 

=  0, 

=  0, 


r= 0 

^e(P)|?=0  =  1' 

The  cathode  scale: 

1  9 

an  +  R  d R 

1  9 


cun  +  R  dR 


{{dll  +  R)ie)  =  3ie(p)|?=1, 

((an  +  R)  =  3  £>J 


r=  1 


(ID 

(12) 

(13) 

(14) 

(15) 

(16) 

(17) 

(18) 


Fig.  1.  Schematic  diagram  of  three  size  scales  within  a  primary  alkaline 
battery  cathode. 


“5  /  £s(p)  9Ce 


96'Mn4+ 


dy 

96*Mn4+ 


‘16 


1 


=  0, 


v=0 


9y 


(3) 


(4) 


es  9r 

0507011 


an  +  $  97? 
1  9 


£>eoo  (an  +  *)  H 


«16  «17  +  ^  9^ 


((an  +  fl)OuCe) 


1 1  Oeoo  8Ce.(P>  -  a5 “7  V(^) Cefp) 


y=i 


‘16 


dr 


r=  1 


(19) 


in  =  “  2 


^Mn44  |  y=i 


1  “3  CMn4+  |  -_j 

1  —  a3 

1  “dC'efp)  \ 


Ce, 


1  —  04  / 


CMn4+l?=0  —  1- 

The  particle  scale: 


(p)exp  (ai4»7(p)  +  «15  (l  -  CW+l^i)) 
exp  (-a14f/(P)  -  ai5  (l  -  CW+j^ 

(6) 


(5) 


dr)  _  an,  /  1 
dR  V^eoo 


i  a 


*2  a-lT^eCp))  -  3 in , 


1  9  9  n 


(7) 

(8) 


+09 


1  + 


“13  |  - 
«ioC 


1  —  a4Ce 


012013017 

“16  (“17  +  R) 
91nae  9Ce 

~dtT~dk’ 


C, 


e  I  R=0 


“5 


9Ce(p) 

9? 


1  9 

r2  dr 


Oeno  r 


9C*e(p) 

dr 


-  3a5a6in 


e\R=0 


0507  9  -□  x 

-^-(r  Ce(p)P(p)), 


(9) 


=  1, 
=  1, 
=  0, 
=  0, 


R=  1 


dfl(p)  -> 

1  o8 

9Ce 

dr 'e(p) 

(cMn4+|.=08)"3 

dR 

=  0, 


R=  1 


+09 


rioCe(p) 


1  04  Ce(p) 


9  In  9e(p)  9Ce(p) 
9Ce(pi  9r 


(10) 


Ce 


?= 0 


=  1. 


(20) 

(21) 

(22) 

(23) 

(24) 

(25) 

(26) 


Hereafter,  we  refer  to  Eqs.  (2)-(26)  as  the  full  cathode  model. 
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Table  1 

Dimensionless  variables 


R-Ri 

Ro-Ri 


^'mh4  i 

CMn4+  — 

LMn4+ 

-  _  Ce(p) 

Ce(p>  _  cj 


Ce 


Ce 

Ce° 

InRi  H  , 


It  - 

DqqO 


I 

D 

D 


i 


eoo 

0~ 

eoo 


e 


?  _  27rR[H(R0  —  R[ )  £emd  . 

hi  —  ;  hi 


tp  - 

1  Vp)  - 


2ttrJj(R0  -  Sj)£EMpf 
ro(VH2o-f°+Ve)/(l-£s(p)) 
2jrRiHF  n 

— T - V 


U(P) 


(  Vh2o  -  *_K+ 


■W1 

2nR\H  (R0  —  Ri )  £emd&2 
7(p)  -  37  ''(Pi 

'  o1 

„  2xRiH(R0  —  Ri)£EMD&2 

"=  - dl - ” 


2jTRiH(R„-  R,)SEMDFClni+ 
»  2nRxH  ( R0  —  Ri)  eemd  . 

Ze(p)  —  —r  n  'e(p) 


I  (l  £s(pj)  r0 


^  _  ft-eoo 

Keoo  —  — fj 

^eoo 


The  relationships  between  dimensionless  variables  in  the 
above  equations  and  their  dimensioned  counterparts  are  listed 
in  Table  1 . 

The  dimensionless  constants  oq-an  appearing  in  Eqs. 
(2)-(26)  are  listed  in  Table  2  along  with  their  typical  value 
ranges.  In  obtaining  the  typical  ranges  for  these  parameters, 
the  discharge  rate  has  been  varied  from  5  to  200  mA  g- 1  of 
EMD,  the  particle  radius  (r0)  from  5  to  250  |xm,  and  the  cath¬ 
ode  thickness  (Ra  —  R,  )  from  0.5  to  5  mm.  Such  variations  are 
indicative  of  the  industrially  relevant  ranges  of  these  param¬ 
eters  [18].  The  remaining  parameter  values  used  to  calculate 
the  typical  ranges  given  in  Table  2  are  given  in  Tables  3  and  4. 
The  dimensionless  constants  used  here  are  similar  to  those 
used  by  Farrell  and  Please  [1]  for  their  particle  scale  model, 
but  here  the  constants  relate  to  the  time  scale  of  cathodic 
discharge  rather  than  the  time  scale  of  particle  discharge. 


Table  2 

Dimensionless  constants  and  their  typical  range  of  values 


Dimensionless  constant 


2nR\H  (R0  -  Ri)£EMDOH+CMn4+F 
2tzR\H(R0  -  R\)  £emeUq 


“3  —  PMn3+CMn4+ 
“4  =  t4Ce° 


2 71-Riff  (Ro  -  RilEEMD^^jFC^ 


Typical  range 


0.032  <  o'!  <1.3 


0.022  <  o2  <  0.87 
1 

0.16 

5.5  x  10-9<<*5< 
0.014 


a6  = 


tU  ru 
K+CMn4+ 

£s(p)  Ce 


0(7  = 


0(  8  = 


(vH2o-i°+t4)  (i-£s(P))c; 


-0 

Mn4+ 


^2  (l  £s(p)) 


£s(p) 


*■0  , 
''•eoo  \ 


ag  = 


A-7T R[H  (Rq  ~  Ri)SEMDRgasTt°+k2 


Frll 


<*10  = 


an  = 


o(  12  = 


<*13  = 


Ce°V h2o 


*K+ 
£EMD  \ 


(l  £s(p; 
*2£EMD 


“14  =  2 


Fril 


“15  = 
“16  = 

“17  = 


AriRiH  (R„  —  Ri)£EMDRgas7?.'2 


2RgaSr 


R0-  Ri 
Ri 

R0  -  Ri 


12 

6.4 

9.7  x  104 

3000  <  “9  <  7.6x  109 
0.74 

0.029 

2800 

8.3  x  10-10 

2.9  x  10'11  < 

“14  5  7.2  x  10~5 

6.8 

2.0  x  10~4  < 

“16  <  0.5 

2000  <  an  < 

2.0  x  104 


3.  Simplified  cathode  model 


We  now  systematically  reduce  the  full  system  of  model 
equations  introduced  above  in  order  to  obtain  a  simplified 
model  that  retains  the  essential  physical,  chemical  and  elec¬ 
trochemical  phenomena  that  comprise  the  discharge  of  a 
primary  alkaline  battery  cathode. 

The  initial  boundary  value  problem  given  by  Eqs.  (2)-(4) 
and  (6)  is  simplified  by  applying  the  method  of  Laplace 
transforms  (See  for  example,  Trim  [19])  to  obtain  a  single 
expression,  which  governs  the  concentration  distribution  of 
Mn4+  within  an  EMD  crystal.  Details  are  given  in  Appendix 
A,  where  we  show  that: 


2  sin  (A ,my) 


C4  (y,  t)  =  1  +  3  /  indt*+TJ2~ 


0  y41-;  sin  (A,„ ) 

m= 1 


X  /  ln 

Jo 


t=t—t*  eXP 


-X.lt* 


at 


d  t*, 


(27) 


Table  3 

Parameter  values  used  in  model 


Parameter 

Values  and  reference 

Cj?  (mol  cm-3) 

0.009  [22] 

CMn4+  (mol  cm-3) 

0.0486  [15] 

Dh+  (cm2  s_1) 

1  x  10-15  [30] 

F(C  mol-1) 

96485.309  [31] 

(A  cm~2) 

5  x  lO"8  [15] 

h  (V) 

0.35 

k2  (S  cm-1) 

1.5  x  102 

h 

4.328 

Rgas  (J  K_1  moP1) 

8.31451  [31] 

T  (K) 

298.15 

>K+ 

0.22  [32] 

14  (cm3  mol-1) 

17.8  [33] 

Vh,o  (cm3  moP1) 

18.07  [33] 

LMn3+  (cm3  moP1) 

20.576 

LMn4+  (cm3  moP1) 

20.576 

ya  (cm) 

2.6  x  10-6  [22] 

£EMD(p) 

0.9  [22] 

“g  (S  cm-1) 

7  x  106  [31] 
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Table  4 

Model  input  data  for  simulation  of  experimental  set-up  used  by  Williams 
[22] 


Parameter 

Value 

Discharge  rate  (mAg-1  of  EMD) 

50 

Particle  radius,  rQ  (cm) 

10  x  10“4,  45  x  10~4, 

100  x  10~4  and  180  x  10“4 

Inner  radius,  R[  (cm) 

1000.0 

Outer  radius,  R0  (cm) 

1000.5 

Height,  H  (cm) 

6.603  x  10“4 

Total  mass  of  cathode  (g) 

4.5 

Mass  of  EMD  in  cathode  (g) 

1.0 

Mass  of  graphite  in  cathode  (g) 

2.9 

where  t*  is  a  dummy  variable  and  the  values  of  Xm  (in  =1,2, 
. . .,  oo)  are  given  by  the  positive  roots  of 

tan(Am)  -  Xm  =  0.  (28) 


Since  all  spatial  variation  of  electrolyte  on  the  particle 
scale  is  eliminated,  Eq.  (10)  becomes: 

=  *e(P)  ( ( CMn4+  |  y— O.g)  +  '  (30) 

Combining  Eq.  (30)  with  the  overpotential  boundary  con¬ 
dition  given  by  Eq.  (15),  we  find  that  at  r  =  0: 

fef  p)  =  0.  (31) 

It  now  becomes  evident  from  the  similarity  of  Eqs.  (7)  and 
(8),  and  Eqs.  (14)  and  (31),  that: 

*e(p)  =  »§)  (32) 

Noting  this,  Eqs.  (17)  and  (18)  now  have  the  same  source 
term  and  they  have  identical  boundary  conditions  at  R  =  1 . 
Thus, 

?e  =  va.  (33) 


Expression  (27)  replaces  Eqs.  (2)-(4),  and  thereby 
accounts  for  diffusion  within  the  crystals  in  the  model. 

An  asymptotic  solution  for  the  Mn4+  concentration  has 
been  obtained  by  Farrell  and  Please  [1],  and  is  based  on  two 
assumptions.  The  first  is  that  a\  is  small.  The  dimensionless 
constant  oq  represents  the  ratio  of  the  time  it  takes  a  pro¬ 
ton  to  diffuse  a  distance  y0  to  the  time  required  to  discharge 
completely  the  cathode.  The  assumption  that  cq  —*■  0  seems 
appropriate  for  all  but  high  discharge  rates  (see  Table  2).  The 
second  assumption  is  that  the  transfer  current,  in,  does  not 
change  on  the  timescale  of  crystal  scale  diffusion.  This  par¬ 
ticular  assumption  is  suited  to  particle  discharge,  because  the 
reaction  distribution  within  a  single  particle  varies  less  when 
the  particle  is  discharged  under  galvanostatic  conditions.  In 
contrast  to  this,  within  a  cathode  being  discharged  galvanos- 
tatically,  each  particle  discharges  with  a  varying  current,  and 
thus  the  transfer  current  is  likely  to  vary  on  the  timescale  of 
crystal  scale  diffusion.  We  note  that  for  the  simplified  case 
of  constant  transfer  current,  Eq.  (27)  reduces  to  that  given  by 
Farrell  and  Please. 

Having  obtained  an  expression  for  the  concentration  of 
Mn4+  ions  in  EMD  crystals,  we  now  turn  our  attention  to  the 
particle  scale.  The  dimensionless  constant  a$  represents  the 
ratio  of  the  time  it  takes  electrolyte  to  diffuse  a  distance  r0  to 
the  time  required  to  discharge  completely  the  cathode.  For  all 
but  very  large  particles  under  very  high  discharge  rates,  015  is 
small  (see  Table  2)  and  we  assume  that  as  -4-  0.  In  addition, 
it  is  noted  that  the  constants  a 5  and  aq  are  0(1)  (see  Table  2). 
Utilizing  these  observations  in  a  regular  perturbation  analysis 
of  Eq.  (9),  it  is  found  that  to  leading  order,  the  spatial  variation 
of  the  electrolyte  concentration  on  the  particle  scale  can  be 
eliminated  from  the  model  (see  Appendix  B),  that  is, 


Ce(p)  (R,  r,  t )  =  Ce  ( R ,  t )  . 


(29) 


A  further  consequence  of  the  loss  of  spatial  variation 
within  the  electrolyte  solution  of  the  porous  particles,  as 
given  by  Eq.  (29),  is  that  the  source  term  within  Eq.  (19)  van¬ 
ishes.  This  term  does,  however,  provide  a  vital  link  within  the 
model,  between  the  reaction  rate  on  the  particle  scale  and  the 
reaction  rate  on  the  cathode  scale.  Thus,  rather  than  apply¬ 
ing  the  leading  order  (i.e.,  O(aj))  expansion  for  Ce(p)  given 
by  Eq.  (29),  in  this  instance  a  version  of  Eq.  (9)  is  obtained 
that  is  accurate  to  0(0:5)  (refer  to  Eq.  (63)  in  Appendix  B). 
Manipulation  of  this  equation  yields  an  expression  for  the 
source  term  in  Eq.  (19)  that  preserves  the  link  between  the 
particle  and  the  cathode  scales  (refer  to  Eq.  (66)  in  Appendix 
B).  Substituting  into  Eq.  (19)  gives: 


3  Ce 

(1  +  «n)  — ^  = 
at 

ajau  3 

~(a17  +  H)  dk 


*16 


0:5(0:17  +  R)  3 R 


Deoo(un  +  R ) 


3 Ce 

3  R 


((ai7  +  R)ieCe)  -  3a6aiiie(p)|?=1- 


(34) 


On  the  cathodic  scale,  the  constant  an  represents  the  ratio 
of  the  conductivity  of  the  KOH  electrolyte  to  the  conductivity 
of  the  graphite  phase.  The  data  in  Table  2  show  that  an  is 
extremely  small,  and  thus  in  the  asymptotic  limit,  a  13  — »•  0, 
Eq.  (20)  can  be  simplified  to  yield: 

aioCe 
1  -  a4.Ce 

Eqs.  (5),  (7),  (17),  (27),  (30),  (34)  and  (35)  with  boundary 
conditions  (12),  (15),  (2 1 )— (23)  and  (25)  and  initial  condi¬ 
tions  (6)  and  (26)  represent  a  simplified  model  of  cathodic 
discharge  for  a  primary  alkaline  battery  system. 


3  In  fle  3  Ce 
3 Ce  ~dk' 


(35) 


3  fj  a\2  » 
fib  =  +  “9 
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4.  Solution  of  simplified  equations 

The  simplified  model  equations  are  solved  numerically 
by  applying  a  finite-volume  method  [20]  in  which  they  are 
discretised  over  the  radial  variable  on  either  the  cathode  or 
particle  scale  and  in  time.  The  reaction  terms  are  fully  lin¬ 
earised,  and  all  other  non-linearities  are  dealt  with  by  apply¬ 
ing  a  fixed-point  iteration  method  (see,  for  example,  [21]). 
The  time  discretisation  can  be  adjusted  by  use  of  a  weighting 
parameter  from  fully  explicit  to  fully  implicit.  The  solution 
is  implemented  in  an  extensive  program  using  MATLAB®. 
The  numerical  algorithm  used  to  solve  the  equation  system 
steps  through  time  and  iteratively  refines  the  solution  of  the 
discretised  system  of  equations  at  each  time  step.  If  the  itera¬ 
tion  process  of  the  solution  does  not  converge,  a  smaller  time 
step  is  chosen  until  it  does.  The  adaptive  nature  of  the  time 
stepping  makes  the  method  robust  over  a  wide  parameter 
range.  The  algorithm  terminates  if  the  required  cell  voltage 
or  finishing  time  is  reached  or  if  the  time  step  is  reduced 
too  much. 

The  simplifications  carried  out  in  the  previous  section 
provide  a  substantial  reduction  in  complexity  of  the  model 
system,  which  now  has  only  one  partial  differential  equation, 
as  compared  with  three  in  the  full  model  system.  The 
ensuing  reduction  in  the  number  of  numerical  calculations 
needed  to  solve  the  system  means  that  the  software  that 
implements  the  numerical  solution  is  easily  run  on  a 
standard  desktop  PC.  For  example,  on  a  desktop  machine 
with  a  2.4  GHz  Pentium®  4  processor  and  512  MB  of 
RAM,  a  typical  discharge  simulation  takes  of  the  order 
of  60. 

4. 1 .  Physical  parameters 

A  list  of  the  parameter  values  used  in  the  model  analy¬ 
sis  that  follows  is  given  in  Table  3.  Many  of  the  values  are 
consistent  with  those  given  by  Farrell  et  al.  [15],  but  the  val¬ 
ues  k\ ,  k.2  and  k 3  were  obtained  by  fitting  the  model  to  the 
experimental  data  of  Williams  [22],  as  discussed  in  Section 

5. 

The  constant  Aq(V)  is  the  coefficient  of  the  so-called 
ion-ion  interaction  term  adopted  by  Farrell  et  al.  [15]  in  the 
full  model.  In  dimensioned  form,  this  term  is  given  by: 

T  =  k\{\  —  CMn4+ 1  j,=  ]  ^Mn3+)-  (36) 

T  appears  in  dimensionless  form  in  Eq.  (5)  of  our  sim¬ 
plified  model.  Such  terms  have  been  used  by  several  authors 
[23-26]  in  order  to  overcome  the  fact  that  manganese  oxide 
reduction  is  not  really  a  single-step  process,  as  defined  by 
reaction  (1),  but  rather  a  combination  of  several  individual 
reduction  processes  [26-29]  for  which  appropriate  kinetic 
data  are  as  yet  unknown. 

The  constants  kj  (S  cm- 1 )  and  £3  are,  respectively,  the  pre¬ 
exponential  and  exponential  coefficients  in  the  manganese 
oxide  conductivity  function  of  the  full  model.  In  dimensioned 


form,  this  function  is  given  by: 

0EMD  =  &2  (l  —  £s(p))  (  +N|  ,  (37) 

\CMn4+/ 

and  appears  in  dimensionless  form  in  Eq.  (30)  of  our  simpli¬ 
fied  model. 


5.  Results  and  discussion 

Williams  [22]  galvanostatically  discharged  a  series  of 
planar  EMD  cathodes  in  KOH  electrolyte  at  a  rate  of 
20,  50  or  100  mAg~'  of  EMD.  Each  cathode  con¬ 
sisted  of  EMD  particles  that  were  taken  from  a  specific 
size  fraction,  namely:  2 r0  <  45  |xm,  77  p,m<  2 r0  <  106  |im, 
150  |i,m  <  2r0  <  300  |xm  and  300  pun  <  2r0  <  500  |xm.  A 
comparison  of  the  output  of  the  simplified  model  with  the 
experimental  data  of  Williams  is  given  in  Fig.  2,  and  Table  4 
lists  the  parameter  values  used  in  the  model  to  simulate  the 
experimental  data. 

In  order  to  simulate  the  planar  geometry  of  the  experimen¬ 
tal  cathodes  using  the  cylindrical  model,  the  inner  and  outer 
radii  of  each  cathode  in  the  model  are  artificially  increased 
whilst  still  maintaining  the  experimentally  observed  cathode 
thickness  and  cathode  volume.  This  procedure  causes  the  cur¬ 
vature  of  the  electrode  to  decrease  and  the  model  equations 
developed  in  the  previous  section  for  cylindrical  geometry 
converge  to  the  corresponding  equations  for  linear  geometry. 

From  the  data  Fig.  2,  it  is  seen  that  the  model  results 
compare  well  with  the  experimental  data.  The  simplified 
model  accurately  captures  the  polarisation  effects  seen  in 
the  experimental  results  with  the  cathodic  discharge  time 
increasing  as  the  EMD  particle  size  is  decreased.  This  is  due 
to  more  uniform  reaction  distributions  within  the  smaller 


Fig.  2.  Comparison  between  experimental  discharge  results  (unshaded  sym¬ 
bols)  and  simplified  model  results  (shaded  symbols)  for  planar  cathodes  with 
configuration  given  in  Table  4.  (O)  *o  =  150-250  pm,  (•)  r0  =  180  pm,  (<) 
r0  =  75-150  pm,  (◄)  r0  -  100  pm,  (0)  r0  =  38.5-53  pm,  (♦)  r0  =  45  pm,  (v) 
r0  =  0-22.5  pm  and  (▼)  rQ  =  10  pm. 
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EMD  particles  that  lead  to  a  greater  utilisation  of  active 
material  within  cathodes  containing  these  particles. 

In  simulating  the  data  of  Williams,  we  have  taken  the  value 
of  CMn 4+  that  appears  in  the  manganese  oxide  conductivity 
function  (given  by  Eq.  (37))  to  be  that  at  80%  of  the  radius 
(i.e.,  0.8yo)  of  agiven  crystal.  This  corresponds  to  the  value  of 
CMn 4+  at  y  =  0.8  appearing  in  Eqs.  (10)  and  (30).  The  effect 
of  choosing  various  positions  within  the  oxide  crystals  at 
which  to  take  the  value  of  CMn 4+  in  order  to  calculate  the  man¬ 
ganese  oxide  conductivity  is  shown  in  Fig.  3.  The  simulations 
are  for  the  cathode  manufactured  by  Williams  that  consists 
of  EMD  particles  in  the  size  fraction  150  |xm  <  2 rQ  <  300  |xm 
(refer  to  Table  4).  The  experimental  discharge  result  is  also 
given  in  Fig.  3. 

The  results  in  Fig.  3  indicate  that  taking  the  value  of  CMn4+ 
at  0.8yo  yields  a  theoretical  discharge  that  corresponds  well 
with  the  experimental  result.  Farrell  et  al.  [15]  and  Farrell  and 
Please  [1]  previously  took  the  concentration  value  appearing 
in  the  conductivity  function  to  be  that  at  the  outer  radius  of  the 
manganese  oxide  crystals.  As  Fig.  3  demonstrates,  however, 
this  approach  appears  to  overestimate  the  resistance  experi¬ 
enced  by  electrons  moving  within  the  oxide  used  by  Williams 
and  leads  to  reduced  discharge  times  in  comparison  with  the 
experimental  data.  To  obtain  more  accurate  predictions  of 
the  conductivity  of  non-uniformly  reduced  EMD,  an  in-depth 
study  into  the  current  paths  and  the  connectivity  on  a  crystal 
scale  would  be  needed. 

In  Section  4,  an  expression  was  obtained  for  the  concen¬ 
tration  distribution  of  Mn4+  within  an  EMD  crystal  (i.e.,  Eq. 
(27))  by  applying  the  method  of  Laplace  transforms.  This 
expression  does  not  depend  on  simplifying  assumptions,  such 
as  a\  — »  0,  and  that  i„  does  not  change  on  the  timescale  of  pro¬ 
ton  diffusion.  Nevertheless,  if  we  are  willing  to  admit  these 
two  assumptions,  then  asymptotic  methods  can  be  applied  to 
the  crystal-scale  proton  diffusion  problem  (i.e.,  Eqs.  (2)-(4) 


Fig.  3.  Comparison  between  experimental  discharge  results  (unshaded  sym¬ 
bols)  and  simplified  model  results  (shaded  symbols)  for  a  planar  cathode  with 
configuration  given  in  Table  4  and  r0  =  75-150  pm  when  various  positions 
within  the  oxide  crystals  are  chosen  at  which  to  calculate  the  manganese 
oxide  conductivity.  Conductivity  calculated  at  (A)  y0,  (•)  f).8y0  and  (■) 

0.6yo- 


and  (6))  in  order  to  obtain  approximate  expressions  for  the 
distribution  of  Mn4+  within  EMD  crystals.  The  analysis  fol¬ 
lows  closely  that  reported  by  Farrell  and  Please  [1]  for  the 
discharge  of  porous  EMD  particles,  and  we  find  that  the  0(1) 
expression  is  given  by: 


CW+  (%  r,  t)  =  CMn 4+  (M)  =  1  +  3  [  h  (f,  t)  dr* 

Jo 


and  the  0(oq )  expression  is  given  by: 


Ovin4+  —  1  +  3  /  in  dr*  4-  a\ 


rS'2 


3 

2  io 


2 sin  (/,„>■) 


/  1  eXP  l 

y  ^r„sm(/.„)  V  a  1 


-kit 


(38) 


(39) 


The  discharge  results  (given  in  terms  of  the  fraction 
of  the  theoretical  capacity  of  the  cathode  that  is  used)  of 
the  simplified  cathode  model  at  various  discharge  rates  are 
presented  in  Fig.  4.  Either  Eq.  (27)  or  (38)  is  used  to  model 
the  distribution  of  Mn4+  within  EMD  crystals.  The  cathode 
configuration  used  here  is  for  a  cylindrical  AA-cell,  details  of 
which  are  given  in  Table  5.  At  low  discharge  rates,  the  use  of 
either  Eq.  (27)  or  Eq.  (38)  within  the  model  yields  very  sim¬ 
ilar  discharge  curves,  however,  as  the  current  is  increased, 
significant  discrepancies  between  the  two  models  are 
observed.  To  understand  why  these  discrepancies  occur,  we 
distributions  of  Mn4+  within  an  EMD  crystal  were  obtained 
at  R  =  Ri  and  r=rQ  as  given  by  Eqs.  (27),  (38)  and  (39)  for  a 
discharge  rates  of  20  mAg-1  of  EMD  (Fig.  5a),  50  mA  g  1 
of  EMD  (Fig.  5b),  or  100  mA  g“ 1  of  EMD  (Fig.  5c).  At  a  low 
discharge  rate,  such  as  that  shown  in  Fig.  5a,  the  distribution 
of  Mn4+  within  an  EMD  crystal  is  essentially  independent  of 
crystal  radius  and  the  simplifying  assumptions  that  constitute 
the  asymptotic  solutions,  namely,  04  — >  0  and  that  i„  does 


Fig.  4.  Discharge  results  for  cylindrical  AA-cell  as  given  by  simplified  cath¬ 
ode  model  when  either  Eq.  (27)  (unshaded  symbols)  or  Eq.  (38)  (shaded 
symbols)  is  used  to  model  the  distribution  of  Mn4+  within  EMD  crys¬ 
tals.  (□/■)  20  mAg-1,  (V/V)  50  mAg-1,  (£)/•)  100  mAg-1  and  «>/♦) 
200  mAg-1  of  EMD. 
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Fig.  5.  Mn4+  concentration  distributions  of  a  crystal  at  R  =  A', .  and  r=ra  as  given  by  Eqs.  (27)  (■),  (38)  (▼)  and  (39)  (•)  for  discharge  rates  of  (a)  20mAg  , 
(b)  50 mAg-1  and  (c)  100  mAg-1  ofEMD. 


not  change  on  the  timescale  of  proton  diffusion,  are  well 
supported  (in  fact  the  a\  value  corresponding  to  Fig.  5a  is 
0.17).  Thus,  the  results  predicted  by  Eqs.  (27),  (38)  and  (39) 
correspond  well  in  this  discharge  regime.  When  the  discharge 
rate  is  increased,  as  in  Fig.  5b  and  c,  the  distribution  of 
Mn4+  within  an  EMD  crystal  becomes  quite  non-uniform 
and  the  discrepancies  between  the  predictions  of  Eq.  (27), 
(38)  and  (39)  become  significant.  Indeed,  at  a  discharge  rate 
of  100 mAg-1  of  EMD,  the  assumption  that  oq  — >  0  can  no 
longer  be  supported  and  the  oq  value  corresponding  to  Fig.  5c 
is  0.86.  In  this  regime,  Eqs.  (38)  and  (39)  become  invalid. 


Table  5 

Model  input  data  for  simulation  of  AA-cell  set-up 


Parameter 

Value 

Discharge  rate  (mA  g-1  of  EMD) 

20,  50,  100  and  200 

Particle  radius,  rQ  (cm) 

100 x  10~4 

Inner  radius,  R[  (cm) 

0.45 

Outer  radius,  R0  (cm) 

0.67 

Height,  H  (cm) 

4.04 

Total  mass  of  cathode  (g) 

10.62 

Mass  of  EMD  in  cathode  (g) 

9.24 

Mass  of  graphite  in  cathode  (g) 

0.8 

6.  Conclusions 

A  simplified  model  for  the  galvanostatic  discharge  of 
primary  alkaline  battery  cathodes  has  been  presented. 
The  simplified  model  was  obtained  by  applying  Laplace 
transform  and  perturbation  methods  to  the  discharge  model 
of  Farrell  et  al.  [15].  In  the  analysis  that  ensued,  it  is 
shown  that  the  three  size  scales  used  by  Farrell  et  al.  to 
describe  the  porous  EMD  cathode  can  be  reduced  to  two 
size  scales  without  the  loss  of  generality.  In  addition,  the 
analysis  demonstrates  that  the  time  taken  for  electrolyte  to 
diffuse  into  a  porous  EMD  particle  is  fast  compared  with 
the  cathodic  discharge  time,  and  that  ohmic  losses  within 
the  graphite  phase  of  the  cathode  can  be  considered  to  be 
negligible.  Furthermore,  the  simplified  model  incorporates 
a  closed  form  expression  for  the  distribution  of  Mn4+  within 
an  EMD  crystal  that  is  not  reliant  on  assumptions  that  may 
break  down  at  high  discharge  rates.  Numerical  solutions 
of  the  simplified  model  equations  are  achieved  via  the  use 
of  a  finite-volume  code  written  in  MATLAB®  that  can 
easily  be  run  on  a  standard  desktop  PC.  The  results  of 
the  code  compare  favourably  with  relevant  experimental 
data. 
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Appendix  A 


Here,  we  consider  the  Laplace  transform  solution  [19]  of 
Eq.  (2)  subject  to  the  boundary  conditions  (3)  and  (4)  and  the 
initial  condition  (6).  Taking  the  Laplace  transform  of  Eqs. 
(2)  to  (4),  treating  in  as  a  function  independent  of  CMn 2+  and 
substituting  the  initial  condition  (6),  we  obtain: 


{  ^Mn4+  }  | 


8y 


—  pL  {  C'mii4+  }  1  ’ 


at  y  =  0, 

dL  {  CMn4+  } 
dy 

and  at  y  —  1 , 

9L  {  CMn4+  } 

dy 


=  0 


=  aqL  {?„  } 


(40) 


(41) 


(42) 


where  p  is  the  transformation  variable  and  L  {  y  }  is  the 
Laplace  transform  of  y.  Using  the  transformation: 


L{eM„U  =  !  +  ^. 

p  y 

we  obtain  the  following  equation  set: 

d2f(y,p ) 


dy2 
/( .0,  p)  =  0 
df(lp) 

dy 


=  at  pf  {%  P) , 


/(l,  P)  =  aiL{/„}  . 


(43) 

(44) 

(45) 

(46) 


Solving  Eq.  (44)  subject  to  the  conditions  (45)  and  (46) 
and  substituting  the  result  into  Eq.  (43),  gives: 

,  1  1 _ opL  {?n  }  sinh  (*JaTpy) _ 

Mn4+  p  y  (V»Tpcosh  (V«T ~P)  ~  sinh  ( ~P)) ' 

(47) 


Applying  the  Convolution  Theorem  of  Laplace  transforms 
to  Eq.  (47)  yields: 

CW+  =  1  +  T  in  |?=?_?*L_1  (§  (p)}  d t*.  (48) 

Jo 

Here,  t*  is  a  dummy  variable  of  integration,  L-1  {§(/?)} 
is  the  inverse  Laplace  transform  of  function  g(p),  where  the 


function  g  ( p)  is  given  by: 


HP)  = 


a  1  sinh  oT[py ) 

y  { v^Ipcosh  (y/uip)  -  sinh  (v^Tp) } 


(49) 


To  obtain  the  inverse  Laplace  transform  in  Eq.  (48),  an 
extension  of  Cauchy’s  integral  formula  [19]  is  used  to  write 
the  inverse  transform  as  an  integral  in  the  complex  plane 
along  the  line  Re(z)  =  y,  where  z  is  the  complex  variable. 
This  is  the  so-called  complex  inversion  integral  of  g(p), 
namely: 


1  ry+iP 

L  1  {g  (p)}  =  —  lim  /  exp  (pt)  g(p)  dp.  (50) 
2m  ft— *0O  Jy_i/3 


Applying  the  Residue  Theorem  [19],  the  right-hand  side 
of  Eq.  (50)  may  be  written  as  the  sum  of  the  residues  of  the 
integrand,  exp  ( pt)  g  ( p ),  at  its  singularities,  which  occur  at 
pm  (m  =  0,  1,  2,  . . .,  00).  It  is  noted  that  po  =  0  and  pm  = 
—Xjn/a\  (m=  1,2,...,  00),  where  the  values  of  Xm  are  given 
by  the  positive  roots  of  Eq.  (28).  The  first  residue  is  3  and  the 
other  residues  are  given  by: 


2  sin  (Xmy)  ( 

- exp  -  ; 

y  sm  (Xm)  \a\  J 


(m  =  1,2, ,  00) . 


(51) 


Using  these  residues  in  Eq.  (48),  we  obtain  an  expres¬ 
sion  for  the  concentration  of  Mn4+  within  an  EMD  crystal, 
namely: 


C4  (y,  t)  =  1  +  3 


sin  (Xmy) 
sin  (Xm) 


x 


exp 


-X 


2  t* 


a  1 


dt*. 


(52) 


Appendix  B 


Here,  we  consider  the  perturbation  analysis  of  Eq.  (9) 
subject  to  the  boundary  conditions  (11)  and  (13)  and  ini¬ 
tial  condition  (16),  in  the  limit  when  0:5  ->  0.  We  begin  by 
expressing  Ce(P)  as  an  asymptotic  expansion  in  powers  of  015, 
namely: 

CVp)  ~  Ce,0  +  tt5Ce,l  +  a^Cet2  +  •  •  •  (53) 

where  Ce.«,  Cc.  1 ,  Ce, 2>  ■  •  ■  afe  independent  of  0:5.  Upon  sub¬ 
stituting  this  expansion  into  Eq.  (9),  the  following  set  of  order 
equations  is  obtained: 


1  d 


9Ce,  0 


=  0. 


At  r  —  1,  Ce,  0  —  Cq. 

9Ce  0 

Atr  =  0,  — —  =0. 

dr 


At?  =  0,  Ce,o  =  L 


(54) 

(55) 

(56) 

(57) 
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0(a5) : 

dce  0  1  3  /.  2ace  A 

di  =r>  dr  dr  J  “  ^ 

~^2Ce,  owg,) 

(58) 

At  r  =  1, 

Ce,  1  =  0. 

(59) 

> 

r-F 

II 

O 

dCe  1 
’  =  0. 
dr 

(60) 

At?  =  0, 

Ce,  1=0. 

(61) 

Solving  Eq.  (54)  subject  to  the  conditions  (55)  to  (57) 
gives: 

Ce, q(R,  r,  t)  =  Ce(R,  ?).  (62) 

Thus,  to  O(ofj),  the  electrolyte  concentration  within  any 
porous  EMD  particle  is  equal  to  the  electrolyte  concentration 
at  the  outer  radius  of  the  particle  as  given  by  Ce(R,  ?). 

Substituting  Eq.  (7)  into  Eq.  (58),  multiplying  by  r 2  and 
rearranging,  yields: 


3 

dr 


'  (7  C’e(p)O^p)) 


a5r 


3C'e(p) 


+  a5a6  — (r2ie(p,). 
3 1  dr 


(63) 


Integrating  Eq.  (63)  with  respect  to  r,  evaluating  the  result 
at  r  =  1  and  multiplying  by  015,  gives: 


,  3Ce,  i  ~  m 

^5  I  ^eoo  ^  ^7Ce,0^(p) 


=  <*5 


dr 

1  3Ce,o 
3~Sf 


f=  l 


+  *e(p)  |f_i 


(64) 


Noting,  from  Eq.  (53)  that 

3Ce,  i  9Ce(p)  2\ 

“5  +  0(0-5), 


dr- 


dr 


(65) 


It  is  seen  that  the  left-hand-side  of  Eq.  (64)  is  equivalent 
to  the  source  term  in  Eq.  (19).  Thus, 


A- no -  a5«7Ce(p)C(j!) 


«5 


dr 

1  3Ce(p) 
3  3 1 


r—  I 


+  a6  *'e(p)  | j 


(66) 
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